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Abstract. 

We model the binary black hole system OJ287 as a spinning primary and a non- 
spinning secondary. It is assumed that the primary has an accretion disk which is 
impacted by the secondary at specific times. These times are identified as major 
outbursts in the light curve of OJ287. This identification allows an exact solution 
of the orbit, with very tight error limits. Nine outbursts from both the historical 
photographic records as well as from recent photometric measurements have been 
used as fixed points of the solution: 1913, 1947, 1957, 1973, 1983, 1984, 1995, 2005 
and 2007 outbursts. This allows the determination of eight parameters of the orbit. 
Most interesting of these are the primary mass of 1.84 ■ 10^°Mq, the secondary mass 
1.46 ■ lO*Af0, major axis precession rate 39°.l per period, and the eccentricity of 
the orbit 0.70. The dimensionless spin parameter is 0.28 ± 0.01 (1 sigma). The last 
parameter will be more tightly constrained in 2015 when the next outburst is due. 
The outburst should begin on 15 December 2015 if the spin value is in the middle 
of this range, on 3 January 2016 if the spin is 0.25, and on 26 November 2015 if 
the spin is 0.31. We have also tested the possibility that the quadrupole term in 
the Post Newtonian equations of motion does not exactly follow Einstein's theory: 
a parameter q is introduced as one of the 8 parameters. Its value is within 30% (1 
sigma) of the Einstein's value q = 1. This supports the no — hair theorem of black 
holes within the achievable precision. We have also measured the loss of orbital 
energy due to gravitational waves. The loss rate is found to agree with Einstein's 
value with the accuracy of 2% (1 sigma). There is a possibility of improving the 
accuracy of both quantities using the exact timing of the outburst on 21 July 2019. 
Because of closeness of OJ287 to the Sun (8 — 12°), the observations would be best 
carried out by a telescope in space. 
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Figure 1. The observations of the brightness of OJ287 from late 1800's until today. 



1. Introduction 



The range of celestial mechanics has been expanding from purely Solar 
system studies to binary stars and planetary systems around other 
stars, and in recent years also to binary pulsars. The most recent appli- 
cation of celestial mechanics is in the binary system of two supermassive 
black holes in the quasar OJ287. The identification of this system as 
a likely binary was made as early as 1982 (Sillanpaa et al. 1988), but 
since the mean period of the system is as long as 12 yr, it has taken a 
quarter of century to find convincing proof that we arc indeed dealing 
with a binary system (Valtonen et al. 2008) . The primary evidence for 
a binary system comes from the optical light curve. By good fortune, 
the quasar OJ287 has been photographed accidentally since 1890's, well 
before its discover in 1968 as an extragalactic object. The light curve 
of over hundred years (Figure 1) shows a pair of outbursts at 11 - 14 
yr intervals and the two brightness peaks are separated by 1 - 3 yrs. 
The system is not strictly periodic, but there is a simple mathemat- 
ical rule which gives all major outbursts of the historical record. To 
define the rule, one assumes that a companion orbits the primary in 
a Keplerian orbit. Then demand that an outburst is produced every 
time the companion passes through a constant phase angle relative to 
the primary, and that also another outburst takes place at the opposite 
phase angle. Due to the nature of Keplerian orbits, this rule cannot be 
written in a closed mathematical form, but its consequences are easily 
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calculated. One of the consequences is that two outburst peaks arise 
per period. By clioosing an optimal value of eccentricity (whicli turns 
out to be 0.7) and by allowing the semimajor axis of the orbital ellipse 
to precess in forward direction at an optimal rate (which turns out to 
be 39° per period), the whole historical outburst record of OJ287 is 
well reproduced. What is more important, the model is able to predict 
future outbursts. The prediction of the latest outburst on 13 September 
2007 (Valtonen 2007,2008) was accurate to one day, leaving little doubt 
about the capability of the model (see Figure 2; the predicted time 
refers in this case to the start of the rapid flux rise, i.e. the beginning 
of the phase 2 of the outburst, as described in the next section). 

The more advanced versions of the model make use of the Post 
Newtonian terms in the two-body orbit calculation as well as likely 
astrophysical details of the radiation mechanism of the outbursts. A 
natural way to induce two outbursts per orbit is to consider an accretion 
disk around the primary body, and to associate the two opposite phase 
angles with impacts of the secondary on this disk. By using the standard 
a disk theory of Shakura and Sunyaev (1973) and its extension to 
magnetic disks by Sakimoto and Coroniti (1981), the problem remains 
mathematically well defined. In fact, using only 6 outbursts as fixed 
points in the orbit, it is possible to solve five orbital parameters in the 
first approximation (Valtonen 2007). The success of this model in pre- 
dicting the 2007 outburst has encouraged us to add the remaining well 
defined outbursts in the model (9 in all), and to solve the problem for 8 
parameters. In this paper we first describe the method of solution, and 
then list the orbital parameters. Finally, we show that the parameters 
fit very well with what we would expect of binary black holes, and show 
that the likelihood of finding an OJ287-like object in the survey of the 
whole sky is of the order of unity, i.e. we have not been particularly 
lucky in the discovery of this binary system in the nearly whole sky 
search of periodic variability of quasars which was initiated by Tuorla 
observatory in 1980. 

As a side result, we may test the idea that the central body is 
actually a black hole. One of the most important characteristics of 
a black hole is that it must satisfy the so called no-hair theorem or 
theorems (Israel 1967, 1968, Carter 1970, Hawking 1971, 1972; see 
Misner, Thorne and Wheeler 1973). A practical test was suggested by 
Thorne and Hartle (1985) and Thorne, Price and Macdonald (1986). In 
this test the quadrupolc moment Q of the spinning body is measured. 
If the spin of the body is S and its mass is M, we determine the value 
of q in 
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For black holes g = 1, for neutron stars and other possible bosonic 
structures q>2 (Wex and Kopeikin 1999, Will 2008). 



2. OJ287 Outburst Structure 



Besides the quasipcriodic pattern of outbursts at about 12 year in- 
tervals, there is another crucial piece of information which helps in 
accurate orbit determination. The radiation at the outburst peaks is 
most likely thermal bremsstrahlung radiation from gas at about 10^ 
K temperature. This is evident from the optical-UV spectral energy 
distribution (Ciprini &: Rizzi 2008) as well as from the unpolarized 
nature of the excess radiation at the outbursts. In contrast, the optical 
radiation of OJ287 is non-thermal highly polarized synchrotron emis- 
sion at 'normal' times, with the optical-UV spectral index of about 
-1.5. The reason for the thermal bremsstrahlung is thought to be an 
impact of a secondary black hole on the accretion disk of the primary 
(Lehto & Valtonen 1996, Sundelius et al. 1997). An outburst begins 
when a bubble of gas torn off the accretion disk expands, cools and 
at some point becomes transparent at optical wavelengths. When such 
a bubble is viewed from a distance, the emission is seen to grow in 
a specific way as the observational front advances into the bubble. 
The size of the bubble, and thus also the rate of development of the 
outburst light curve, is a known function of distance from the center 
of the accretion disk. The process has three stages: (1) a spherical 
bubble becomes transparent at its forward side which causes the initial 
rise in the light curve. (2) When the bubble is seen approximately 
half-way through, it is fully transparent, and the flux rises sharply in 
the light-crossing time of the bubble. (3) The fully transparent bubble 
continues its expansion adiabatically, and consequently the radiation 
falls in power law manner. Figure 2 illustrates the three stages for the 
observations of the 2007 outburst, with the dash-dot and dashed lines 
outlining the first and the third stages of the outburst, respectively. 
From the point of view of timing, it is important to identify the initial 
moment of transparency when stage (1) begins. For the 2007 outburst 
this moment was at late hours of 10 September GMT; expressed as 
fractional years the moment is 2007.692 (the line drawn earlier than 
this has no physical significance in Figure 2). Similar timing fits have 
been carried for all 9 outbursts, and the corresponding timings are 
listed in Table 1, including the measurement uncertainties. Note that 
the timing ranges are slightly narrower than in Valtonen et al. (2010) 
where the spin determination was first carried out. It makes the finding 
of solutions harder but not impossible. 
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Figure 2. Observations of the 2007 outburst (crosses) together with the hnes mod- 
eUing the phase (1) and phase (3) evolution of the radiating bubble. Phase (2) 
evolution happens on 13 September when the brightness increases sharply. 



3. PN-accurate orbital description 

Wc calculate the binary orbit using PN-accurate orbital dynamics that 
includes the 3PN- accurate conservative non-spinning contributions, the 
leading order general relativistic, classical spin-orbit and radiation re- 
action effects (Barker and O'Connell 1975, Damour 1982, Kidder 1995, 
Mora and Will 2004). The relevant PN-accurate equations of motion 
can be written schematically as 

X = = XO + XlPiV + X50 + XQ 

+X2PAr + X2.5PAr + X3P7V , (2) 

where x = xi — X2 stands for the center-of-mass relative separation 
vector between the black holes with masses mi and m2 and xq repre- 
sents the Newtonian acceleration given by xq = — x; m = mi + m2 
and r = |x|. The PN contributions occurring at the conservative IPN, 
2PN, 3PN and the reactive 2.5PN orders, denoted by xipjv, X2Piv , 
xapjv and X2.5PAr respectively, are non-spin by nature, while X50 is the 
spin-orbit term of the order 1.5PN. 

The employed PN-accurate equations of motion are in harmonic 
coordinates and satify covariant spin supplementary condition. In the 
present analysis, we also included the precessional motion for the spin 
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Table I. Outburst 
times with estimated 
uncertainties. The 

1982/3 outburst is used 
as a reference time. 



1912.970 


± 0.010 


1947.282 


± 0.0005 


1957.080 


± 0.020 


1972.940 


± 0.005 


1982.964 


± 0.0 


1984.130 


± 0.002 


1995.841 


± 0.0005 


2005.740 


± 0.005 


2007.692 


± 0.0005 



of the primary black hole due to leading order general relativistic spin- 
orbit coupling. 

The quadrupole-monopole interaction term xq, entering at the 2PN 
order, reads 



XQ =qX 2c'^A 



5(n • si) 



n-2(n -81)81 >, 



where the parameter q, whose value is 1 in general relativity, is intro- 
duced to test the black hole 'no-hair' theorem. The Kerr parameter x 
and the unit vector 81 define the spin of the primary black hole by the 
relation Si = Gmfxsi/c and x is allowed to take values between 
and 1 in general relativity. The unit vector n lies along the direction of 

X. 

The initial orientation of the binary orbit is perpendicular to the 
accretion disk of the primary. The spin axis of the primary black hole 
is tilted 6 degrees relative to the spin axis of the accretion disk, and 
remains at this angle while the black hole spin axis precesses around the 
disk spin axis, with a period of 1300 yrs. The disk itself also wobbles, 
but with a smaller amplitude and with a period of 120 yrs. This shows 
up as modulation in the long term evolution of the light curve. For 
details, see Valtonen et al. (2009, 2010). 



4. Timing experiments 

We now require that the correct orbit must reproduce the nine outburst 
times of Table 1 within the stated error limits. An automatic search 
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algorithm has been constructed which after many iterations finds a 
solution, starting from given initial values. The algorithm performs 
orbit integrations over and over again until the correct orbit is found. 
Because of the tolerance in the observed outburst times, the parameters 
of the solutions also have their associated ranges. Table 2 lists these 
parameters and their ranges. The latter are the three sigma error limits 
of the orbital parameters. 

Tabic II. Solution parameters. 





39°. 1 ± 0°.l 


mi 


(1.84 ±0.01) • 10^" Mo 


m2 


(1.46 ± 0.1) • 10® Mo 


X 


0.28 ±0.03 


(t>o 


56°.3 ±1°.0 


eo 


0.658 ± 0.001 


q 


1.0 ±0.9 


td 


0.74 ± 0.04 



We will now discuss the physical significance of some of the param- 
eters. The precession rate 39°. 1 per period is defined as the average 
change of the apocenter phase angle over the last 200 yrs. During each 
orbit the time and the phase angle is noted down when the two bodies 
are at their maximum separation. The exact values are determined by 
extrapolation between integration steps. The value given in Table 2 is 
the average increase of the phase angle during the 200 yr of integration, 
starting in 1856. A good first estimate is given by the timing of the 2005 
outburst which came a year earlier than expected in non-precessing 
models. The precession rate leads to the value of the primary mass 
mi = 1.84 ■ 1O^°M0 which is quite consistent with the values of quasar 
black hole masses in a large sample of SDSS quasars; the masses extend 
well beyond 2 • IO^^Mq (Vestergaard et al. 2008). 

As to the expected frequency of such high mass values in a mag- 
nitude limited quasar sample, we may make the following estimate. 
OJ287 is highly variable mostly due to strong beaming. At its faintest 
it goes to niB = 18 which may be taken as its intrinsic (unbeamed) 
brightness. There are about 2 • 10*^ quasars in the sky brighter than 
this magnitude limit (Arp 1981). Many of these quasars host binary 
black holes, perhaps as many as 50% (Comerford et al. 2009). Thus 
potentially there are about 10^ bright binary quasars in the sky to be 
discovered. Vestergaard and Kelly (private communication) estimate 
that there is a 25% probability that about one hundred of them are 
more massive than 10^° Mq. It is not unlikely that one of them would 
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possess a jet which points more a less directly toward us, and that 
therefore the quasar would appear as a blazar. The brightness of the 
quasar is magnified by orders of magnitude by the jet, and so is the 
probability of its detection. 

On the other hand, most of the binary black hole systems should 
have longer periods than OJ287 which is in the last stages of inspiral, 
with only 10^ yr left out of its potential 10'' yr lifetime (Volonteri et 
al. 2009). This compensates at least partly for the increased discovery 
rate due to beaming. Thus the number of bright short period binaries in 
the 10"*^*^ Mq category which could be discovered by techniques similar 
to the discovery of OJ287 is of the order of unity. We may have been 
lucky to discover one such example, OJ287, but anyhow this is only 
an order of magnitude calculation, and more OJ287's may turn up in 
future searches of periodically varying quasars. 

What about the secondary mass? Its value was calculated by Lehto 
and Valtonen (1996) clS Tft'2 1.4 • 10^ M0, updated to today's Hubble 
constant and to the outburst peak brightness observed in 2007. This 
estimate is based on the astrophysics of the disk impacts and on the 
amount of radiation produced in these events, and thus it is totally 
independent of the orbit model. This agrees well with our more exact 
value m2 = 1.46 • 10*^ Mq. The fact that the mass ratio is as high as 126 
is quite natural since equal mass binaries do not easily merge within the 
Hubble time (Makino & Funato 2004) and also because the accretion 
disk would be unstable if the mass ratio were less than 100. 

The initial apocenter eccentricity at the beginning of our simulation 
is Co = 0.66. It refers to the osculating Keplerian orbit and the oscu- 
lating eccentricity value varies considerably over the orbit. If defined 
by using pericenter and apocenter distances as for a Keplerian orbit, 
the eccentricity is e = 0.7. The initial eccentricity at the beginning of 
the final inspiral, say, 10^ yr prior to merger, must have been e 0.9, 
a reasonable value at the initial stage of merging black holes (Aarseth 
2008) . Then by today the eccentricity would have evolved to its current 
value. This is typical behaviour of binary black holes of large mass ratios 
in galactic nuclei. 

The parameter gives the ratio of the viscosity parameter Ug to 
the accretion rate in Eddington units. The value td = 0.74 implies 
14 lb 0.5 for this ratio. The mass accretion rate may be about 0.005 of 
the Eddington rate (Bassani et al. 1983) which gives ag = 0.07 ± 0.03, 
a reasonable value. 

The dimensionless spin parameter is x = 0.28. Some of the outburst 
timings are particularly sensitive to x, among them the 1973 outburst. 
Figure 3 illustrates how the outburst light curves depend on x in this 
case. Thus we infer that the primary black hole spins approximately 



IAU261ValtonerL.tex; 8/01/2010; 14:23; p. 8 



Black hole spin 



9 



28 




1972.94 1972.96 1972.98 1973 1973.02 1973.04 1973.06 1973.08 1973. 

Year 



Figure 3. The brightness measurements of OJ287 in 1972/3 (dots with error bars). 
The error bars rejjresent short-time variability, to lesser extent measurement errors. 
The three lines show theoretical light curves for three values of x- 

at one quarter of the maximum spin rate allowed in general relativ- 
ity. There is an additional observation which supports this spin value. 
OJ287 has a basic 46ib3 day periodicity (Wu ct al. 2006), which may 
be related to the innermost stable orbit in the accretion disk. However, 
since we presumably observe the accretion disk almost face on, and 
there is an m = 2 mode wave disturbance in the disk, this is likely 
to refer to one half of the period. Considering also the rcdshift of the 
system, and the primary mass value given in Table 2, this corresponds 
to the spin of x = 0.35 ±0.06 (McClintock et al. 2006). For comparison, 
it has been estimated that x 0.5 for the black hole in the Galactic 
center (Genzel et al. 2003). We will have a unique opportunity to im- 
prove the spin determination in 2015. If the spin is x = 0.28 the next 
outburst begins on 15 December, 2015. The range of possible outburst 
times extends from 26 November if % = 0.31 to 3 January 2016 if 
X = 0.25. 

The values of q cluster around q = 1.0 with a standard deviation of 
0.3. Figure 4 illustrates the 1995 outburst and how its timing depends 
on q. The error limits of q may be further narrowed down in 2019, to 
the level of 20% (1 sigma) if the timing is carried out accurately. The 
fact that q = I with the current accuracy indicates that the 'no-hair' 
theorem is valid. It is interesting that our result converges at the proper 
value for a general relativistic black hole. 
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Beginning of 1995 outburst 
at 1995.841 



q=1.3,--;. 



q=0.7 



1 1.5 
Day - Juiian Day 2450026.86 



Figure 4- Observations of brightness of OJ287 in 1995 (dots with error bars). The 
three hnes are theoretical curves for the rising part of the flux (for q = 0.7, 1.0, 1.3). 



We have also tested the sensitivity of the solution to the 2.5PN 
(radiation reaction) term. One may multiply the 2.5PN term by 1.0 it 
0.06 and still find a solution. Thus we confirm general relativity with 
6% (3 sigma) accuracy. It is possible to improve this accuracy also by a 
factor of 2 in the timing of the 2019 outburst. The main observational 
problem with this outburst is that it occurs close to the conjunction of 
OJ287 with the Sun. The separation in the sky is only 12° at the start 
of the outburst on 21 July 2019, and decreases to 8° during the week of 
the main burst. It is an observational challenge, probably too much of 
a challenge for ground based observations. Thus a space born telescope 
measurement would be desirable. This is a unique opportunity to test 
the 'no-hair' theorem and the radiation reaction term which will not 
be repeated for many decades. 
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